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Abstract 

In principle, the state space of a chaotic attractor can be partially or wholly 
reconstructed from interspike intervals recorded from experiment. Under cer¬ 
tain conditions, the quality of a partial reconstruction, as measured by the 
spike train prediction error, can be increased by adding noise to the spike 
creation process. This phenomenon for chaotic systems is an analogue of 
stochastic resonance. 
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The seminal articles |l| and the theorem of Takens in ||^ demonstrated the theoretical and 
practical possibility of reconstructing the topological structure of the state space underlying 
an experimental system, using the measurement of a generic scalar or multivariate signal 
from the system. This possibility is especially welcome for nonlinear systems, where the 
potential exists for extremely complicated state space attractors. A great deal of subsequent 
research effort has gone into developing data processing techniques for the detection, analysis 
and exploitation of nonlinear, and in particular chaotic, processes. 

Often, a delay coordinate reconstruction of a compact attractor from an evenly-sampled 
time series of measurements can be found that is topologically equivalent to the attractor. 
The function mapping the attractor to the reconstructed copy is called an embedding. Ac¬ 
cording to 0, as long as the embedding dimension is greater than twice the box-counting 
dimension of the attractor, an embedding results for a probability-one choice of measurement 
functions. Recent theoretical work has attempted to widen the scope of dynamical data that 
can lead to an embedding. It has been shown, for example, that a similar reconstruction 
result holds when using spike train data (the recorded times between firings) from a model 
integrate-and-fire dynamical system with chaotic dynamics [Q. 

Attractor reconstruction can be viewed as a type of information transfer. In the case 
of a topological embedding, no information is lost. The set of states of the underlying ex¬ 
perimental system is reproduced exactly in the copy that is reconstructed from measured 
data. In other cases, the reconstruction may be incomplete. When a chaotic signal is fed 
into a threshold crossing detector, the time-delay plot of time intervals between crossings 
reconstructs something akin to a Poincare section of the underlying chaotic attractor. The 
dimension is decreased by one 0. A recent study 0 of neuron models subjected to chaotic 
input points out that the two-dimensional FitzHugh-Nagumo differential equation (here¬ 
after referred to as FHN2) acts as a kind of threshold-crossing “hlter” for input signals, and 
in particular fails to completely reconstruct the attractor which generated the input signal. 
One focus of the present article is to clarify the distinction, for attractor reconstruction 
purposes, between threshold-crossing (TC) and integrate-and-hre (IF) hlters. In contrast to 
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the fact that FHN2 acts as a TC hlter, we exhibit a variation on FHN2 which acts princi¬ 
pally as an IF hlter instead, and does fully reconstruct the input attractor. This differential 
equation model, which we shall denote FHN3, is a three-dimensional version of excitable 
FitzHugh-Nagumo-type dynamics. With appropriate parameter settings, it has the follow¬ 
ing remarkable property: When a signal from a system attractor is added to one of FHNS’s 
variables, another variable undergoes a deterministic spiking behavior whose interspike in¬ 
tervals, embedded as m-tuples, embed the original attractor in m-dimensional reconstruction 
space. Just as Takens’ theorem guarantees that multidimensional state space information 
can be condensed into a single evenly-spaced time series, the same can be accomplished with 
spike timings from the FHN3 hlter. 

Viewing attractor reconstruction as information transfer raises questions about the effi¬ 
ciency of the transfer process, and the possible ehects of noise on this process. Surprisingly, 
in some instances noise can have a benehcial ehect on information transmission, analogous 
to the stochastic resonance phenomenon |p] observed for multistable potentials and excitable 
media. This seemingly contradictory ehect is the observed amplihcation of a hltered signal 
achieved when stochastic noise is added to the input signal. 

Stochastic resonance has been shown to amplify signals generated by linear models, such 
as sine waves. Recently, it has been shown P] that the same basic ehect can be seen with 
aperiodic (stochastic) input signals if the means of measurement is appropriately modihed. 
Since the spectrum of a random signal is not discrete, the SNR must be replaced in this case 
with a power norm sensitive to shape-matching and/or signal correlation. 

Our present interest in stochastic resonance is somewhat diherent. Our goal is in max¬ 
imizing the amount of state information carried by the spike train. We pose the question 
whether the quality of state information from a deterministic signal (as opposed to the 
stochastic realization used in |P) can be improved by injecting random noise into the pro¬ 
cess. In particular, we are interested in the case where the input signal is chaotic. To 
determine the quality of state information contained in the spike train, we measure the 
ability to predict the spike train from its own history, using nonlinear prediction techniques. 
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We will show examples in which the interspike interval prediction error of the output signal 
decreases (predictability increases) with increasing noise added to the input signal. Showing 
that nonlinear predictability is enhanced by adding noise is analogous to the enhanced SNR 
shown in stochastic resonance studies. 

The hlters we will use to create spike trains capable of carrying low-dimensional deter¬ 
ministic state information are based on the well-known FitzHugh-Nagumo equation |0. The 
two dimensional system FHN2 

ev = —v{y — 0.5)(n — 1) — w + S 

w = V — w — h (1) 


is a simple differential equation that exhibits a fast spike followed by a refractory period. 
There is an equilibrium at {v,w) = {vo,Vo — b), where Vq is a real-valued root of — 
0.5)(vo — 1) = S + b — vq. a stability check of the equilibrium vq shows the existence 
of a supercritical Hopf bifurcation for Sh = vh{vh — 0.5)(u/f — 1) + vh — b, where vh = 
0.5 — \/3 — 12e/Q. Therefore, if 6, e are hxed and the bifurcation parameter S is increased, 
the system undergoes a Hopf bifurcation at Sh, resulting in a periodic orbit of the system 
encircling the formerly stable equilibrium. The periodic orbit is manifested in rhythmical 
spiking by the variable v. For example, setting 6 = 0.15, e = 0.005, there is a Hopf bifurcation 
point at Sh ~ 0.112331.... For S < Sh, the system is quiescent; the equilibrium is stable. 
For S > Sh, the system spikes at a rate of approximately 1 Hz. 

Now consider FHN2 as a nonlinear hlter by substituting for the constant S' in (]^) a signal 
S(t) from another system. Fig. 1 shows a plot of the variable v from (]l|) where S has been 


replaced by a signal from the Rossler system |]l0 


X = T{-y - z) 
y = T(x + ay) 

z = r(6 + {x — c)z) (2) 


where the standard parameters are set to a = 0.36, b = 0.4, c = 4.5, and r = 0.5 causes 
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the trajectory to run at half speed. The input signal, also plotted in Fig. 1, is S{t) = 
0.09 + 0.013a;(t), where x{t) is the x variable of (^. The bias of the signal is {S{t)) = 0.093, 
and its root mean square amplitude is {{S{t) — {S{t))y) = 0.035. Fig. 2 demonstrates the 
threshold-crossing detection capability of FHN2. When the peak height of S{t) is greater 
than 0.15088, FHN2 hres a burst of spikes. Note that this threshold is signihcantly higher 
than the Hopf bifurcation value Sh, which would be the threshold in the limit of an S{t) 
which oscillates inhnitely slowly. 

Although we expect the spike sequences to carry state information of the Rossler system, 
because of its threshold detection behavior we do not expect it to carry enough to reconstruct 
the entire attractor. On the other hand, the beneht of a TC hlter is that noise can in some 
cases enhance the reconstruction quality, as measured by prediction error. As in studies of 
stochastic resonance, we will add white noise to the input signal of the hlter (in this case, 
the FHN2 spike generator). The equation with noise term is 

ei) = —v{v — 0.5)(u — 1) — w + S{t) + 

w = V — w — h, (3) 


where e = 0.005, b = 0.15, and ^{t) is Gaussian white noise with zero mean and autocor¬ 
relation (^(t)^(s)) = 2DS{t — s). For small values of the noise level D (including all those 
considered here), the variable v exhibits a clearly distinguishable spiking behavior, often in 
bursts of more than one spike, as in Fig. 1. For analysis purposes, we found it more con¬ 
venient to collect series of interburst intervals, each dehned to be the elapsed time between 
the hnal spike of one burst and the hrst spike of the next burst. After using (^ to make 


a series of 1024 interburst intervals , we used a standard nonlinear prediction algorithm 
to measure the level of determinism in the series. The fact that state information from a 
deterministic system is contained in a spike train, even when the spike train is chaotic, can 
be detected by measuring the nonlinear predictability of the interburst intervals. If it can 
be shown that the ISI series is predictable “beyond the power spectrum”, that is, if there 
is predictability beyond that which is guaranteed by linear autocorrelation, then there is 
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evidence of nonlinear dynamics in the series. 


The prediction algorithm works as follows. Given an ISI vector Vq = (tjg,... ,tio-m+i)) 
the 1% of other reconstructed vectors 14 that are nearest to Vq are collected, omitting vectors 
14 close in time. The ISI for some number h of steps ahead are averaged for all k to make 
a prediction. That is, the average is used to approximate the future interval 

tig+h- The difference p — ti^+h is the h-step prediction error at step 4- We could instead 
use the series mean m to predict at each step; this h-step prediction error is m — tjo+fc. 
The ratio of the root mean square errors of the two possibilities (the nonlinear prediction 
algorithm and the constant prediction of the mean) gives the normalized prediction error 
NPE = {{pio—tio+hYy^'^/{{'n^ — tio+hYy^^ where the averages are taken over the entire series. 
The normalized prediction error is a measure of the (out-of-sample) predictibility of the ISI 
series. A value of NPE less than 1 means that there is linear or nonlinear predictability in 
the series beyond the baseline prediction of the series mean. 

The results of the predictability of the interburst interval series from (|^) are shown in 
Fig. 3. For these parameter settings, unlike those for Fig. 1, no spikes occurs in the absence 
of noise. As the noise power D is increased from zero, spikes begin to occur for very small 
noise levels, although the interburst series show no predictability (NPE 1) until D is raised 
beyond 10“^^. The prediction error then drops to a minimum and raises again when the 
noise becomes large enough to swamp the system. The clearly noticeable improvement in 
predictability due to extremely small noise input is essentially a stochastic resonance effect. 
These results show evidence of nonlinear determinism, since Gaussian-scaled surrogate series 
created from all burst series considered in Fig. 3 have NPE 1. 

A slight alteration in the FitzHugh-Nagumo equations yields a nonlinear hlter that acts 


as an integrate-and-hre processor. Dehne the system FHN3 by 


ii = —au — cw + S{t) 
ev = —v{v — 0.5) (u — 1) + u — dw 

w = — w — b. (4) 
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This system is similar to FHN2 in that if S{t) is set to be a constant parameter S, there 
is a Hopf bifurcation as S is increased. Setting parameters a = 0.1,6 = 0.15, c = 0.5, d = 
0.5, e = 0.005, the bifurcation point is Sh ~ —0.059. Its success as an information processor 
is shown in Fig. 4. As with FHN2, we replace the parameter S with an input signal S{t) 
from the Rossler attractor. The signal is S(t) = 0.0023a;(t) — 0.04, which corresponds to bias 
{S(t)) = —0.04 and root mean square amplitude \J{{S(t) — {S(t))y) = 0.006. Comparing 
with FHN2 in Fig. 1, we see a marked difference in the way FHN3 processes the input signal. 
Fig. 5(a) shows a three-dimensional plot of the vectors (tj, fi_i, tj_ 2 ), where ti = Ti — Ti_i 
is the time interval between spikes of the v variable of FI1N3. Fig. 5(b) shows a similar 
reconstruction where the input signal S{t) is the x-coordinate from the Lorenz equations 
1^ i: = T{a{y—x)),y = T{px—y—xz), z = T{—j3z+xz), where the parameters are set to the 


standard values a = 10, p = 28,/3 = 8/3, and r = 0.01. Apparently, the interspike intervals 
recovered from (^) do an effective job of reconstructing the chaotic attractor which produced 
the input signal S{t), for both the Rossler and Lorenz examples. Nonlinear prediction on 
a length 1024 series of spikes created as in Fig. 4 yields NPE = 0.1. This very low NPE 
supports the visual indication in Fig. 5(a) of a faithful reconstruction of the underlying 
Rossler attractor. This is similar to the mechanism that was studied in the generic integrate- 
and-fire model of @], where firing times Tj were generated recursively by 


i-i + 1 


S{t)dt = 0 


( 5 ) 


for a hxed threshold 0. Theoretical reconstruction results for spike trains generated by 
model d^) are discussed in Q. 

Creating spikes using FI1N2 or FI1N3 means imposing a type of highly nonlinear hlter on 
the attractor signal, a hlter which edits out amplitude information (since the spike waveforms 
are essentially alike) and converts the information entirely to event timings. Our purpose 
is to gain insight into the data processing methods used in systems which communicate 
through spike timings, as is conjectured for certain neural systems |^. We have shown by 
example that noise may be useful for this communication, in that it can amplify transmission 
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of deterministic, nonlinear state information as measured by nonlinear prediction error. For 
the latter spike generation model (FHN3), we have the possibility of complete reconstruction 
of attractor states. 
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FIGURES 


FIG. 1. The solid curve is the variable v of FHN2, the FitzHugh-Nagumo equation (|I|) with S 
replaced by S(t) = 0.09 + 0.013x(t), where x(t} is a solution of the Rossler system (^). The dashed 
curve is S{t). When a peak of S(t) is greater than « 0.15, a burst is triggered in FHN2. 

FIG. 2. Peak heights of the signal S{t) from Fig. 1 graphed versus time. The height is plotted 
as an asterisk if it triggers a burst from FHN2; as an open circle if not. All of the asterisks he 
above all of the open circles, signifying precise threshold detection by FHN2. 

FIG. 3. Normalized prediction error of spike trains generated by (|^, where S{t) is a signal 
formed using the Rossler x-variable of (^) with bias 0.075 and rms amplitude 0.020 (open cir¬ 
cles) or 0.023 (asterisks). As the input noise power D increases, the NPE displays a minimum, 
corresponding to maximum information transfer. Each plotted point is an average over 5 noise 
realizations; standard error is less than 0.02 for each. 

EIG. 4. The solid curve is the variable v of EHN3, equation (Q), with S(t) = 0.0023a:(t) — 0.04. 
The dashed curve is x(t), the x-variable of the Rossler attractor (^). A plot of 3-tuples of interspike 
intervals from this equation is shown in Fig. 5a. 

FIG. 5. Interspike interval reconstructions of (a) the Rossler-FHN3 intervals from Fig. 4 (b) 
Lorenz-FHN3 intervals from (|^ with S(t) = .0005x(t) — 0.04, where x(t) is the x-variable of the 
Lorenz system. In (b), fewer points are plotted, and they are connected with line segments. 
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